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Abstract 

We examine the axial U{1) symmetry near and above the finite temperature phase transition 
in two-flavor QCD using lattice QCD simulations. Although the axial U{1) symmetry is always 
violated by quantization, i.e. the chiral anomaly, the correlation functions may manifest effective 
restoration of the symmetry in the high temperature phase. We explicitly study this possibility by 
calculating the meson correlators as well as the Dirac operator spectral density near the critical 
point. Our numerical simulations are performed on a 16^ x 8 lattice with two flavors of dynamical 
quarks represented by the overlap fermion formalism. Chiral symmetry and its violation due to the 
axial anomaly is manifestly realized with this formulation, which is a prerequisite for the study of 
the effective restoration of the axial U{1) symmetry. In order to avoid discontinuity in the gauge 
configuration space, which occurs for the exactly chiral lattice fermions, the simulation is confined 
in a fixed topological sector. It induces finite volume effect, which is well described by a formula 
based on the Fourier transform from the ^-vacua. We confirm this formula at finite temperature 
by calculating the topological susceptibility in the quenched theory. Our two flavor simulations 
show degeneracy of the meson correlators and a gap in the Dirac operator spectral density, which 
implies that the axial f/(l) symmetry is effectively restored in the chirally symmetric phase. 

PACS numbers: 

Keywords: Lattice QCD, Topology, overlap-Dirac operator, U(l), symmetry 



I. INTRODUCTION 

Chiral symmetry plays a key role in understanding the nature of the vacuum of Quantum 
Chromodynamics (QCD). With Nj flavors of massless quarks, the QCD lagrangian has a 
symmetry U{Nf)i ® U{Nf)Ft under field rotations involving flavor mixings for left-handed 
and right-handed quarks independently. Among them the flavor- singlet part of the chiral 
(or axial) transformation, which forms U{1)a, is violated by quantization. This is known 
as the chiral or axial anomaly. At zero temperature, the remaining symmetry SU{Nf)L ® 
SU{Nf)ji (g) U{l)v is spontaneously broken down to SU{Nf)Y (g) f/(l)y , and Nj — 1 massless 
Nambu-Goldstone bosons appear. At finite temperature, in the massless limit, a restoration 
of symmetry is expected back to SU^Nf)^ ® SU^Nf)^ (g) U{l)v above a critical temperature 
T 

■^ c- 

This symmetry breaking and its restoration are characterized by a vacuum expectation 
value of a flavor-singlet scalar operator qq, which is called the chiral condensate —{qq) = S, 
the order parameter of this phase transition. According to the Banks-Casher relation [l| 

S = 7rp(0), (1) 

non-zero chiral condensate implies an accumulation of low-lying eigenvalues of the Dirac 
operator D. Here, p(A) is a spectral density of D defined as p(A) = (l/V) J2x'{^i^ ~ ^')) 
with V a four-volume of space-time and (■ ■ ■ ) denotes an expectation value. The Banks- 
Casher relation ([1]) is satisfied in the thermodynamical limit, i.e. infinite volume limit, 
followed by the limit of vanishing quark mass. This implies that the density of low-lying 
eigenvalues of D must disappear above the transition temperature. 

The chiral anomaly is reflected in the particle spectrum (in the low temperature phase) 
as a mass-splitting of the corresponding pseudo-scalar meson, the t] meson in the case of 
Nf = 2 for example, from the flavor non-singlet pseudo-scalar particles, i.e. n mesons. In 
terms of (valence) quark flow diagrams, the difference between flavor- singlet and non-singlet 
mesons comes from a disconnected diagram in the meson two-point correlation functions. 
It is found that the main contribution to the disconnected diagram is from low-lying quark 
eigenmodes if one decomposes the quark propagator into the contributions of individual 
eigenmodes of the Dirac operator (2|. An interesting question then emerges: above the 
finite temperature phase transition, where we expect the suppression (or even absence) of 
the low-lying eigenmodes, does rj become degenerate with flavor non-singlet pseudo-scalars? 



If so, the axial U{1) symmetry is effectively restored at least in the particle spectrum at 
high temperature. This may have significant impact in the nature of the finite temperature 
phase transition, since the symmetry controls the order and critical exponents of the phase 
transition |3|. ^ 

Such argument was originally made by Cohen ^], followed by several theoretical studies 
9|. Numerical study using lattice QCD, on the other hand, has been missing, because 
most numerical simulations at finite temperature performed so far used the staggered fermion 
formulation, which realizes chiral symmetry only partially. Even more important problem of 
the staggered fermion is in the so-called rooting procedure, i.e. a square-root of the fermion 
determinant is taken to represent two flavors of light quarks in the nature out of four artificial 
flavors built in the staggered fermion. The fiavor-singlet axial transformation is ill-defined 
in such a theory. There is a recent study of the Dirac spectral density using the staggered 
fermion [10|. In this paper, we present the first numerical study of this problem employing 



the overlap fermion formulation [UJ, [12] , for which the full flavor and chiral symmetries are 
realized on the lattice including the flavor-singlet axial symmetry existing at the Lagrangian 
level and violated by quantization. 

Although the numerical cost for the simulation of the overlap fermion is substantially 
higher than that of the staggered fermion, its simulation has become feasible by recent 
developments of machines and algorithms, and in fact has been extensively performed by 
the JLQCD and TWQCD collaborations in the past years [iSj. The overlap fermion is 
constructed with the overlap-Dirac operator Dqv, 

Dov = mo[l + 75Sgn(i7iy(-mo))], (2) 

where H]v{—fno) is the hermitian Wilson-Dirac operator Hw^—mo) = 'y^D^^^—rrio) with 
a large negative mass —mo, which is of order of the lattice cutoff. This operator satisfies 
the Ginsparg- Wilson relation [ij], -DqvTs + 75-Dnv == -Dov75-Dov/"^0; through which the exact 
chiral symmetry on the lattice can be defined 15|. The operator ([2]) has a singularity due 
to the sign function, when an eigenvalue of H^y^—niQ) crosses zero as the background gauge 
field varies. This corresponds to a boundary between two adjacent topological sectors of 
gauge field. Numerical simulation with currently available algorithms becomes prohibitively 
costly when one crosses this singularity; the JLQCD and TWQCD collaborations took a 
strategy to fix the topological sector during the Monte Carlo simulation. Physical quantities 



at the ^ = vacuum is reconstructed later by correcting the finite volume effect of 0(1/1^) 



due to fixing the topology [l6|, Il7 |. 



Since the artifact due to fixing topology is essentially a finite volume effect, the same 
strategy should work at finite temperature as far as the spatial volume is sufficiently large. 
Detailed discussions are given in Section [ITl We numerically check this property by calcu- 
lating the topological susceptibility on quenched lattices, which is given in Section lllli As 
discussed in |17|, even when the global topology is fixed, local topological fiuctuations are 
still active and the topological susceptibility can be extracted from long-range correlation 
of topological charge densities. The results at finite temperature obtained in this work are 
compared with conventional calculations from the fiuctuation of global topology found in 
the literature. 

We then discuss the results of our exploratory simulations of two-fiavor QCD at finite 
temperature using the overlap fermion formulation in Section IIVI The global topological 
charge is fixed to zero. We analyze the low-lying eigenvalue spectral density of the overlap- 
Dirac operator, as well as the meson correlators to investigate the effective restoration of the 
axial f/(l) symmetry. Both connected and disconnected parts of the correlator are recon- 
structed from the low-lying eigenmodes, which are calculated and stored in advance. The 
disconnected part produces the difference between the meson channels related by the f/yi(l) 
symmetry transformations. We show that, in the high temperature phase, the disconnected 
contribution is indeed vanishing towards the chiral limit of sea quarks. 

A similar study of the restoration of the axial U{1) symmetry has recently been made by 
the HotQCD collaboration |l18|] using the domain-wall fermion formalism. Like the overlap 
fermion, the domain-wall fermion realizes exact chiral symmetry on the lattice, but only in 
the limit of infinitely large fifth dimension. In practical simulations, the chiral symmetry 
is slightly violated and an additive mass renormalization rrires appears. The size of Tribes 
is typically an order of MeV, and poses a significant problem when one tries to identify 
the near-zero eigenvalues, which are of the same order or even lower. The overlap fermion 
employed in this work enables a clear identification of the near-zero modes and their effects 
on the physical quantities. 

Reports of our work at earlier stages are found in 



191,120 



II. PHYSICS AT FIXED TOPOLOGY 

In QCD, the fluctuation of global topology Q is necessary to guarantee the cluster de- 
composition property, which is one of the fundamental conditions necessary in constructing 
meaningful quantum field theory. It is easy to see that the different topological sectors have 
';o be added with a weight of the form e*^*^, through which the QCD ^-parameter is defined 
2l|. The formulae to relate the physical quantities obtained in a fixed topological sector 
Q to those in the = vacuum are developed in [l6|, [iTJ. They are valid also at finite 
temperature as outlined below. 

Let Z{9) be a partition function of QCD with the 0-term put in a volume V^ = L^ x /5. 
Here, L and /3 are spatial and temporal extent of the box respectively (/3 = Nto). The 
partition function Zq at a fixed global topological charge Q can be written as a Fourier 
transform of the 0-vacua: 

Zq = ^ I deZ{e)exp{teQ) = ^ [ de exp{-l3F{e)L^), (3) 

where F{e) = Eo{e) - i9Q/V. 

Near zero temperature, Eq{6) represents the vacuum energy density, Z{0) = 
exp{—f3EQ{9)L^). The integral in (|3]) may be evaluated using the saddle point expansion 
around 

e. = ^^{l + 0{6% (4) 

where Xt is defined through an expansion of Eq{6) in terms of 6, 

oo 



and S = Q/ixtV)- The result up to small 6"^ terms is given by 



1 f Q' 

:exp 



wxr W 



l-:.^ + 0(^,<52 



(6) 



The global topological charge distributes as Gaussian with a width specified by the topo- 
logical susceptibility xt- At finite volume, the correction appears at order 1/V. The basic 
assumption in obtaining ([6]) through the saddle point expansion is that Eq{9) takes its 
minimum value a.t 6 = and is analytic around there. 

At non-zero temperature, /3 plays the role of inverse temperature. As the first excited 
state energy £'1(6') becomes significant compared to the temperature, (3{Ei{9) — Eq[9))L^ ~ 



1, the partition function Z{6) receives contributions from the excited states as Z{6) = 
'^j^exp{—f3En{0)L^), and the relation ([3]) has to be modified to include them. This can be 
done by redefining the "energy" as f3E{6)L^ = —\nZ{6) and accordingly the "free energy" 
/3F{9)L^ = f3E{9)L''^ - i9Q/V. As far as the analyticity property of E{9) around 6* = is 
unchanged, the same saddle-point approximation can be applied, and (jS]) is valid with Xt 
and C4 evaluated at finite temperature. The condition that l/V"^ and 6"^ have to be small in 
([6]) is still applied. Since /3 is fixed, the condition implies large spatial volume L^ to make 
V = (3L^ sufficiently large. 

The same conclusion is obtained by considering the "energy" defined through a transfer 
matrix in one of three spatial directions. In this way, the spatial extent L plays the role 
of inverse temperature and the whole argument given at vanishing temperature remains 
unchanged. Again, the volume V appearing in the formulae is the space-time volume /3L^. 

Along the same line of argument, one can obtain the relation between correlation functions 
at fixed Q and those at fixed 9. For instance, for a correlator G{9) that is CP-even at 6' = 0, 
one obtains llj 



Gq = G{0) + G^^\0) ^ 



2xtV 



Q' C4 



XtV 2xlV 



+G^'\^)^.+0{V-'). (7) 

Namely, the correlator at a fixed topology Gq is written in terms of G'(O) and its second 
derivative G^'^\Q) = dG{9)/d9\Q=f), and so on. The first correction is again of the order of 
l/V. 

One important example of ([7]) is that for a two-point correlation function of the topological 
charge density u{x). It can also be written in terms of fiavor-singlet pseudo-scalar density 
operators mP{x) using the fiavor-singlet axial- Ward- Takahashi identity. Fixing the global 
topology, a constant correlation remains at long distances: 

hm (mP(a:)mP(0))Q = U^ - Xt - ^) + 0{e-^\^\). (8) 

i'|->-large V \ V ^Xt>^ J 

At finite temperature, the long distance |x| must be taken in the spatial direction. 

The constant correlation appearing in ([H]) is understood as the effect of fixing topology. 
Let us consider a sector of Q = 0, as an example, i.e. the global topology is constrained 
to zero. If there is a positive topological charge fiuctuation locally near the origin, there 
would be more chance to find negative fiuctuation apart from there in order that (mPlx)) 



is summed up to zero when integrated over space-time. The correlation {'mP{x)mP{0)) is 
thus negative at long distances and is proportional to the ability of having local topological 
fluctuations, which is characterized by Xt- Since the effect must vanish on large enough 
volume, one expects a contribution of the form —Xt/V, 

The relation ([H]) suggests a possibility to extract Xt from the measurement done in the 
flxed topology. Indeed, it was successfully performed at zero temperature in J2|. 

The last term in ([8]) represents a physical correlation due to the flavor-singlet pseudo- 
scalar particle, here denoted as rj. In the lattice calculation of the flavor-singlet correlation 
functions, there appear connected and disconnected quark-flow diagrams. Both diagrams 
contain a slow decay mode due to the flavor non-singlet pion, which cancels between the 
two diagrams, and only the rapidly decaying channel of heavier rj particle remains. In the 
quenched theory, there is no pion component in the disconnected sector and a modifled 
contribution from the so-called hairpin diagram has to be considered, as discussed later in 
Section [mCl 

When Xt is very small, which is expected for Nf = 2 QCD above Tc, our method for 
extracting xt is no longer valid. However, the disconnected propagator is still useful since 
its absence (or presence) itself is a signal for the Ua{1) restoration (or breaking). If there 
is a sizable gap in the Dirac spectral density, it is likely that the disconnected diagram is 
highly suppressed and Xt = 0. Conversely, if the eigenvalue density is non-zero but very tiny, 
Xt is likely to have a non-zero value. In this way, even when xt is small, we can indirectly 
investigate whether Xi = or not, and ?7yi(l) is restored or not, through the flxed-topology 
simulations. 

III. STUDY OF FIXED TOPOLOGY IN QUENCHED QCD 

Before performing two-flavor QCD simulations at flnite temperature, we carry out a 
quenched study in order to validate the strategy of extracting the physics of the 9 = vacuum 
from the flxed topology simulations at flnite temperature. We measure the topological 
susceptibility at flnite temperature using the method outlined in the previous section and 
compare the result with those in the literature obtained with the conventional method of 
calculating the variance of the global topological charge (Q^). 



A. Setup and datasets 

In order to fix the global topological charge throughout the Hybrid Monte Carlo (HMC) 
simulation, the JLQCD collaboration introduced two extra species of unphysical Wilson 

n 

fermions with a large negative mass —mo [22j. With the mass of order of the lattice cutoff, 
the extra degrees of freedom are irrelevant in the continuum limit. Since the overlap-Dirac 
operator D^^ is built upon the hermitian Wilson-Dirac operator Hy/{—inQ) as in ([2]), a change 
of the index of Dqv accompanies a zero-crossing of one of the eigenvalues of Hw{—niQ). The 
extra Wilson fermions generate a fermion determinant of the form dei[Hy/{—mQY] and 
such zero-crossing is prohibited. We also introduce twisted-mass ghosts to cancel unwanted 
contributions from the extra Wilson fermions. The net effect for the Boltzman weight in the 
path integral is 

det[i7vK(-mo)2 + /i2]' ^ ^ 

where /i is the twisted mass given to the ghosts. The eigenmodes of \Hw{—mQ)\ above /i do 
not contribute to the Boltzman weight effectively, and only the near- zero modes are affected. 
In this work, we chose /i = 0.2 in the lattice unit. 

The suppression of the near- zero modes of iJiy(— ttiq) does not spoil local topological 
fluctuations that give rise to the topological susceptibility Xt- Indeed, using the same for- 
mulation, Xt was successfully calculated and conflrmed to be consistent with the expectation 
from chiral perturbation theory |2j. 

We generated flnite-temperature data in the pure gauge theory on a lattice of size 24^ x 6. 
With Nt = 6, the transition temperature corresponds to the lattice spacing a ~ 0.11 fm, 
which is in the region where the locality of the overlap-Dirac operator is satisfled with our 
choice of lattice actions. We take a range of lattice spacing 0.09-0.13 fm, which corresponds 
to /3 = 2.35-2.55 with the Iwasaki gauge action. The lattice spacing is estimated from the 
heavy quark potential measured on independent zero temperature lattices with an input of 
vq = 0.49 fm. For each parameter, we accumulated 0(300-500) configurations separated by 
50-100 HMC trajectories in the trivial topological sector Q = 0. The lattice parameters are 
summarized in Table HI 

The autocorrelation time of standard thermodynamical observables like internal energy 
and the order parameters, the Polyakov loop and chiral condensate, are negligible. Decor- 
relation of topology-related quantities, such as the lowest eigenvalue of the overlap-Dirac 



Type 


/3 


a (fm) 


T (MeV) T/Tc 


-'•eigenval 


-' 'correlators 


FT 


2.35 


0.132 


249 


0.86 


106 


106 


FT 


2.40 


0.123 


268 


0.93 


336 


336 


FT 


2.43 


0.117 


281 


0.97 


101 


101 


FT 


2.445 


0.114 


288 


1.00 


424 


420 


FT 


2.45 


0.113 


290 


1.01 


584 


584 


FT 


2.46 


0.111 


295 


1.02 


251 


245 


FT 


2.48 


0.107 


306 


1.06 


420 


321 


FT 


2.50 


0.104 


316 


1.10 


379 


218 


FT 


2.55 


0.094 


348 


1.20 


487 


487 


CT 


2.58 Q = 


0.099 


331 


1.10 


299 


235 


CT 


2.58 Q = 


1 0.099 


331 


1.10 


262 


257 



TABLE I: Parameters for the pure gauge simulations on a 24^ x 6 lattice with the topology-fixing 
extra Wilson fermions (labeled as FT). The global topological charge is fixed to Q = in these 
runs. Lattice spacing is estimated from the heavy quark potential. The number of configurations 
used for eigenvalue calculations the (disconnected) correlator calculations is given in the columns 
of A'cigenvai and A'^correiators- Also listed are the runs without fixing topology (CT) at /? = 2.58 (with 
the Iwasaki gauge action). The critical temperature for this case is ~300 MeV [23]. The total 
number configurations for this run is 1,069. Out of them, the sectors of global topological charge 
IQI = and 1 are selected for measurements of the disconnected correlators. 



operator, is very 



slow at higher temperatures, which is a well known effect on fine lattices 



(see for example |2J, |25| |) 



We estimate the finite temperature phase transition point T^ by an inflection point of 
the Polyakov loop, which is more precise than its susceptibility with our statistics. The 
corresponding (3 value is /3 = 2.445, which corresponds to 288 MeV. For each lattice, a value 
of T/Tc is also listed in Table [H The precise value for T^ is irrelevant for the conclusions of 
this section. 



10 



B. Eigenvalue spectral density 

Before analyzing the topological susceptibility, we investigate the eigenvalue distribution 
of the overlap-Dirac operator, which is closely related to the topological susceptibility as we 
discuss below. Above the transition temperature the density of near- zero eigenvalues must 
vanish as dictated by the Banks-Casher relation jl|^. On the other hand, a sharp peak at 
very low eigenvalues has been observed in previous works at or above the critical temperature 
[Sl, |26|, l27[. These were conjectured to be related to a presence of local topological object. 



such as an instanton-anti-instanton pair. If it is the case, these near-zero eigenmodes should 
have significant contributions to the topological susceptibility. 

Eigenvalues of the overlap-Dirac operator ([2]) are calculated applying the Implicitly 
Restarted Lanczos algorithm for a chirally projected operator P^DovP+, where P+ = 
(1 + 75)72. It gives a real part of the eigenvalue of D^v, ReAov, and then the complex 
Aov is reconstructed by solving a relation |1 — Xov/mo\'^ = 1, which is a direct consequence of 
the Ginsparg- Wilson relation. The complex eigenvalue Aov is projected onto the imaginary 
axis as A = ImAov/(l — ReAov/2mo). 

When applying the overlap operator ([2]) , we approximate the sign function using a rational 
approximation with the Zolotarev coefficients after subtracting a few lowest-lying eigenmodes 
of the kernel operator Hw{—mQ). Taking the degree of the rational function to be 16th, the 
precision of the approximation is kept better than 10~^°. With this tiny error we confirmed 
that the low-lying eigenvalues of Dov are obtained to 8-digits or better. The precision of the 



approximation is kept better than 10~^. Further details may be found in 13|. 

In Figured] we plot the eigenvalue spectral density observed on pure gauge lattices. In this 
study, S is defined with the lattice regularization and not converted to the continuum regu- 
larization schemes, such as the MS scheme. On the quenched lattices, the thermodynamical 
limit is equivalent to the infinite volume limit, and our lattices are already sufficiently close 
to this limit when extracting S from the data. (Note that there is no long-range mode in 
the pure gauge theory.) 

From the plot, we find that at low temperature {/3 = 2.35-2.43) the eigenvalues accumulate 



^ Strictly speaking, this is valid only when the chiral condensate is an order parameter of the phase 
transition. In the quenched theory, this is not absolutely necessary, but we assume it as a working 
hypothesis in this section. 
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FIG. 1: Eigenvalue spectral density in quenched QCD at finite temperature. The histograms 
plotted with filled symbols are those below the critical temperature. Those above the critical 
temperature are plotted with open symbols. The data at would-be critical point (/? = 2.445) 
are shown without symbols. The error bars associated with the symbols are estimated using the 
jackknife method. The inset shows a magnification of the near-zero mode region, showing an 
accumulation of very low eigenmodes. 

to the amount corresponding to vrp(O) ~ (200-300 MeV)^. There is no sharp change observed 
at the critical temperature {/3 = 2.445), but p(0) decreases towards higher temperature and 
eventually vanishes at /3 = 2.55. 

More interestingly, we confirm the presence of the peak near A = at around the critical 
temperature. Note that our configurations are generated in the Q = sector and there is no 
exact zero modes. The near- zero modes responsible for the peak have very small eigenvalues 
(less than 20 MeV) but they are still non-zero. We will show that these near-zero modes 
give the main contribution to the topological susceptibility. 
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C. Disconnected correlation functions 



So far, very few groups have studied the topological susceptibility Xt at finite temperature 



in the pure SU(3) gauge theory. The most complete results are found in two papers |28l. |29|. 
Both calculated the variance of the global topological charge to obtain Xt = {Q'^)/V- The 
former used the geometrical definition of the topological charge, i.e. a discretized version 
of F^yF^y, on configurations generated with the standard Wilson action. The topological 
susceptibility was shown to be stable around the zero-temperature value (~ (180 MeV)^) 
until the transition temperature, where it starts decreasing. A disadvantage of the geomet- 
rical method for measuring the topological charge is that it requires some cooling steps that 



potentially affect the final results. The other work 29| used the Liischer-Weisz action and 
directly counted the number of zero modes of a lattice Dirac operator, which approximately 
satisfies the Ginsparg- Wilson relation to measure the global topological charge through the 
index theorem. The results of these two works reasonably agree with each other. We take 
the numbers from 29|] for a comparison with our results. 

We calculate the topological susceptibility on the gauge configurations generated at a fixed 
global topological charge using the formula ([H]) • It requires a calculation of the disconnected 
quark-fiow diagram, appearing in the evaluation of the correlation function {mP{x)mP{0)). 
Namely, we have to calculate 

D,,{x, 0) = (Tr [-f,S{x, x)] Tr [75^(0, 0)]) (10) 

for each space-time point x. Here, S{x, y) is a quark propagator obtained by an inversion 
of the overlap-Dirac operator Dqv for a given source point y. Since the numerical cost for 
a direct calculation is too expensive, one typically uses some stochastic techniques. In this 
work, on the other hand, we introduce a representation of S{x, y) in terms of the eigenvalues 
Afc and eigenvectors ipk{x): 

k '^ 

This is an exact representation when the sum is taken over all eigenmodes. We truncate the 

sum to include only the low-lying modes, which should dominate the long-distance correla- 
tions. We validate this approximation by inspecting the correlation functions constructed 
with 30, 40 and 50 low-lying eigenpairs at /3 = 2.40. Each eigenpair includes two terms in 
( fTTj) : one from the eigenmode of A^ and ipk{x), and the other is that of A^ and 'J5'il'k{x), 
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FIG. 2: Disconnected correlator approximated by different numbers of eigenmodes, i.e. 30, 40 
and 50 eigenmodes. Data at /3 = 2.40 and m = 0.05, which is below the transition temperature. 
The continuous line is a fit of the 50 eigemodes data. 

which is also an eigenmode of Dov For the valence quark mass m in f lTTj) we take am = 0.01, 
but the final result for the disconnected diagram is almost independent of am. 

Figure [2] shows the disconnected correlator with different numbers of low-lying eigenpairs. 
We find that the disconnected correlation function is indeed dominated by the low-lying 
eigenmodes and there is no significant difference between 40 and 50. Therefore we can 
safely assume that the evaluation of the disconnected quark diagram using the 50 eigenpairs 
is sufficiently precise and adopt this procedure in the following analysis. 

On the quenched lattices, we only evaluate the disconnected quark-line contributions to 
the correlator, since the long-distance correlation from the pion channel does not exist in 
the disconnected contribution. Instead, there is a so-called "hairpin" contribution [30| of 
the form 



fP ~^~, 2""^0 -^-, 2 fP ' 



(12) 



where mT^ is the pion mass and fp denotes the matrix element to annihilate a pion to 
the vacuum through the pseudo-scalar density operator. The singlet mass-parameter mo 
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represents a coupling between quark loops in the quenched vacuum. In the coordinate 
space, it corresponds to the functional form ~ /|,mo(l + mTrt) exp(— -m^t) for zero spatial 
momentum. We use this function to fit the lattice data at finite temperature, taking t in the 
spatial direction, together with a constant term —Xt/V representing the fixed topology effect. 
We neglect the subleading effect of —C4/{2xtV'^) in (jS]). The pion mass m^ is effectively 
determined by combining a fit of the connected diagram with that of the disconnected 
diagram. The overall coefficient, such as fpfriQ is treated as a free parameter. 

An example of the fit is shown in Figure [21 which is for a lattice slightly below the 
critical temperature. It demonstrates that the ansatz describes the lattice data well. The 
underlying assumption of our analysis is that the pion channel gives dominant contribution 
at long distances also at finite temperature. This is a reasonable assumption below the 
critical temperature Tc, and seems to be valid even above Tc as our data are well fitted. 

D. Topological susceptibility 

The results for the topological susceptibility Xt obtained at a fixed global topology Q = 
from the constant long-distance correlation of the flavor- singlet correlator are shown in 
figure O Our data are plotted by black dots, which are in good agreement with those from 



291] (squared and diamonds) obtained from the Q^ distribution. 

In order to further cross-check, we also accumulated 1,069 configurations without fixing 
topology by eliminating the extra Wilson fermions iQ. This run is carried out slightly above 
the transition temperature by choosing the lattice spacing a ~ 0.10 fm. In Table [H it is 
denoted as CT. The result for xt obtained by counting the number of exact zero modes is 
plotted in Figure |3] by a cross at around T/Tc ~ 1.1, which shows a good agreement with 
our determination from the fixed topology run. 

It is also interesting to see the consistency of this topology-changing run by selecting 
configurations of a given Q and analyzing them with the method for a fixed topology. We 
pick up two subsets of configurations with global topological charge \Q\ being and 1, and 
calculate the disconnected correlation function to extract Xt- The results are plotted in 
Figure [3] at around T/Tc ~ 1.1 with different symbols. They show perfect agreement with 
the standard method, i.e. extracted from {Q'^)/V, as well as with the result of the fixed- 
topology run. This provides a firm numerical evidence that the method to extract xt at 
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FIG. 3: Topological susceptibility in quenched QCD calculated on a 24^ x 6 lattice. Data points 
(black dots) are obtained from a fit of the correlation functions (see text). The errors are statistical 



only (jacknife binned). Reference data points from [29|] (but on 20'^ x 6 and 16'^ x 6 lattices) are 
shown by diamonds and squares. The data coming from selected sectors in a changing topology 
run are also shown (Q = 1 is slightly shifted for readability). 

topology-fixed configurations works as expected, at least for Xt > 100 MeV. 

We also check that the eigenvalue spectral density does not depend on whether or not 
the configurations are generated at a fixed topological sector. The eigenvalue distribution 
is compared in Figure IH At a matched lattice spacing a ~ 0.1 fm, thus at a matched 
temperature T/T^ :i^ 1.1, the spectral density at Q = is essentially unchanged even when 
the configuration is generated with the constraint for the global topological charge. 

IV. DYNAMICAL QCD WITH OVERLAP FERMIONS 

Given the theoretical formulation and numerical validations for the strategy to extract 
the 6 = vacuum physics from the fixed-topology simulations, we embark on a dynamical 
simulations of finite temperature QCD using the overlap fermion formulation. The lattice 
size in this exploratory study is 16^ x 8 and the global topological charge is fixed to Q = 0. 
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FIG. 4: Comparison of the eigenvalue spectral density from the runs with and without fixing 
topology at a matched lattice spacing a ~ 0.1 fm. The fixed-topology run at /3 = 2.50 (black, 
thin line) is one of those shown in Figure [TJ The result of the Q = configurations selected out 
of the run without fixing topology is overlayed (red, dashed line). These two results are in good 
agreement despite the slight mismatch of the lattice spacing (thus temperature). The spectral 
density of all configurations in the topology-changing run (thick line) show a slight deviation at 
the lowest bin. This is understood as an eff'ect of the exact zero modes on the Q ^ configurations, 
that are taken out in this plot. They repel the nearby eigenvalues and make the spectral density 
lower in the vicinity of A = 0. 

A. Run parameters 



In Table [TTl we list simulation parameters. We use the Iwasaki gauge action together 
with the extra Wilson fermions and associated ghosts ([9]). The temporal size A'"t = 8 is 
chosen to generate smooth enough gauge configurations at around the critical temperature 
to guarantee the locality properties of the overlap operator 3l|- The aspect ratio L/Nt is not 
sufficiently large for finite temperature simulation. This is a possible source of systematic 
errors especially in the vicinity of the phase transition, but we do not consider such errors in 
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13 am a (fm) T (MeV) T/T^ iVcigcnval A'correlators 

2.18 0.01 0.144 172 0.95 118 100 

2.18 0.05 0.144 172 0.95 350 320 

2.20 0.01 0.139 177 0.985 187 187 

2.20 0.025 0.139 177 0.985 303 272 

2.20 0.05 0.139 177 0.985 279 279 

2.25 0.01 0.128 192 1.06 335 331 

2.30 0.01 0.118 208 1.15 512 479 

2.30 0.025 0.118 208 1.15 226 183 

2.30 0.05 0.118 208 1.15 281 281 

2.40 0.01 0.101 243 1.35 477 319 

2.40 0.05 0.101 243 1.35 210 210 

2.45 0.05 0.094 262 1.45 80 

TABLE II: Parameters for finite temperature QCD simulations with two flavors of dynamical 
overlap fermions. The lattice size is 16'^ x 8. The global topology is fixed to Q = 0. 

this work, which is the first attempt to extract the physics related to the axial f/(l) sector 
from the overlap fermion simulations. 

Our lattice volume and statistics are not sufficient to precisely determine the transition 
point solely from the generated configurations. Instead, we assume a value of the critical 



temperature T^ = 180 MeV in two-f 
existing analysis at zero temperature 



avor QCD. The lattice spacing is estimated from an 



32| . From its /3 value (/3 = 2.30), the lattice spacing in 
the range (/3 = 2.18-2.45) is obtained assuming the renormalization group running. System- 
atic uncertainty associated with this estimate exists, but it does not affect the conclusion of 
this paper, which is qualitative. 

The quark mass for degenerate up and down quarks is taken in the range 0.01-0.05 in 
the lattice unit. This corresponds to the bare mass of 14-70 MeV for the lattice close to the 
transition point (/3 = 2.20). 
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B. Dirac operator spectral density 

As in our quenched analysis, we calculate the low-lying eigenvalues of the overlap-Dirac 
operator on two-flavor QCD ensembles. Since the chiral condensate {qq) gives an order 
parameter of the finite temperature phase transition, according to the Banks-Casher relation 
(II]) the spectral density near the zero eigenvalue provides a direct measure of the phase of 
the system. However, we should note that the relation holds only in the thermodynamical 
limit (V — )■ oo then iriq — )■ 0), and the results at finite volume and quark mass have to be 
taken with care. 

For our main interest in this paper, i.e. the effective restoration of the axial f/(l) 
symmetry, the spectral density plays an unique role. Let us consider the susceptibilities 
Xs = J d^x{Jsi^)Jsi^)) ^^^ Xtt = J d'^^{Jni^)J'^i^)) of iso-triplet scalar and pseudo-scalar 
operators Js{x) = q{x)T"'q{x) and j"(x) = q{x)'j5T°'q{x). (r" is the Pauli matrix to specify 
the isospin component. We use a somewhat old notation '5' for the iso-triplet scalar state. 
In the modern terminology, it is called oq.) Using the property that the eigenmodes of the 
Dirac operator appears as a complex-conjugate pair of A^ and A^ and with their eigenvectors 
simply related by 75, i. e. ipk and 'j5ipk, one can show that the difference of the susceptibihties 
is written in terms of the eigenvalues: 

Xn-Xs= / dXp{X) (13) 







[m 



2 + A2)2- 



Disappearance of Xtt — Xs suggests the effective restoration of the axial f/(l) above the 
critical temperature, or at least the anomalous violation of the axial f/(l) cannot be seen in 
this iso-triplet (pseudo-) scalar channel. 

In the broken phase, i.e. p(0) 7^ 0, the difference [T3] diverges as 1/m in the chiral limit, 
which is understood as a contribution from the long-distance correlation due to the Nambu- 
Goldstone pion channel. In the symmetric phase p(0) = 0, on the other hand, the difference 



survives under 
is not fulfilled 



;he condition p(A) ~ A" with a < 1. It was recently shown that this condition 

9[, i.e. a> 2. 

For other channels, the axial f/(l) restoration cannot be simply parameterized by only 
using the spectral function and the details of the eigenvectors are relevant. We could still 
expect the important role played by the near-zero modes, and it is important to identify the 
strength of its suppression, i.e. the power a. There is even a possibility to find a gap in 
p(A), i.e. zero density from A = up to some value Ac. 
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For this reason, we study the eigenvalue spectral density also in two-flavor QCD. Since 
it requires the inflnite volume limit followed by the chiral limit, more investigation would 
be necessary for conclusive results. The flrst results obtained in this work with exact chiral 
symmetry would still give valuable information towards this goal. 

On the conflgurations generated as listed in Table [ITl we calculate 50 lowest eigenvalues 
and associated eigenvectors of the hermitian operator ■y^Dov Paired eigenmodes of Dqv 
can be reconstructed from them, and we effectively have 100 low-lying eigenmodes. The 
numerical method is the same as the one employed in the quenched study. Since the global 
topological charge is flxed to zero, we do not have exact zero-modes. 

We plot the results in Figure E] after rescaling them to a dimensionful unit to compare 
the spectra at different temperature values. We investigated the Monte Carlo histories of 
the lowest mode and found no evidence to suggest long auto-correlations. Also, doubling 
the statistics for a couple of temperatures did not signiflcantly change the initial result. 

First of all, at low temperature (/3 = 2.18, T ~ 170 MeV, top panel of Figure [5]) we 
flnd signiflcant number of near- zero eigenmodes suggesting S ^ 250 MeV as known in the 



zero temperature case 
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34| . Such near-zero modes found at a heavier sea quark mass 
am = 0.05 are suppressed at am = 0.01. This is what should happen on a finite volume 
lattice, even though we expect non-zero density at A = in the infinite volume limit. Such 



finite V and finite m^ scalings are studied in detail at zero temperature 33N36|. 

Near the transition temperature (/3 = 2.20 and 2.25, T ~ 180-190 MeV, middle panel 
of Figure [5]), the result is qualitatively unchanged. Even above the transition temperature 
(/3 = 2.30 and 2.40, T ~ 210 MeV, bottom panel of Figure [5]), we still see a similar number 
of near-zero modes when the quark mass is large {am, = 0.05). This is consistent with our 
observation in the quenched theory (Figured]); it indicates that at the quark mass am = 0.05 
the system is qualitatively similar to the pure gauge theory. 

Once the quark mass is decreased towards chiral limit at higher temperatures, not only the 
near-zero modes, say those below 10 MeV, but the modes up to ~40 MeV disappears, which 
is very different from what we found on the lower-temperature lattices (T ~ 170 MeV). Our 
observation suggests the following picture: In the thermodynamical limit, p(0) disappears at 
the critical point (by definition). At a temperature slightly above, the suppression towards 
the chiral limit occurs more rapidly and at some higher temperature a gap opens. This means 
a stronger suppression than any power a of the form ~ A". Unfortunately, any quantitative 
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FIG. 5: Spectral density of the massless overlap-Dirac operator in two-flavor QCD. Top and 
bottom panels are the data clearly below and above the critical temperature, respectively. The 
middle panel corresponds to those around the transition point. The jackknife errors are shown for 
each bin of the histogram. When the histogram is terminated at the lower end, it implies that we 
find no eigenmode below that value. The statistical error in that case is also zero, because we use 
the jackknife method. The lighter the color the lighter the mass. 

argument about the power a and the point where gap opens would not be possible with the 
currently available data. There is even a possibility that the gap develops right above the 
critical point. Much more extensive data at several quark masses and volumes would be 
necessary for a definite conclusion on this point. 
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C. Meson correlators 

Further information about the restoration of the symmetry can be extracted by directly 
inspecting the degeneracy of meson correlators. Under the flavor- singlet (or iso-singlet) axial 
f/(l) transformation, the pion channel is related to 6, likewise the a channel is related to 77. 
On the other hand, flavor non-singlet (or iso-triplet) chiral transformation connects a to vr 
and 1] to 6. Therefore, above the transition temperature we expect a degeneracy between 
a and vr as well as between rj and S. If the axial U{1) symmetry is effectively restored, 
we should see the degeneracy between a and rj as well as between vr and 6. Namely, all 
four-channels should become degenerate. 

Since the difference between vr and rj or between a and 6 comes from disconnected quark- 
flow diagram, the U{1) restoration means the absence of the disconnected diagram. The 
disconnected contributions to the iso-singlet scalar and pseudo-scalar correlators are written 
as 

Dix,y) = {TT[S{x,x)]TT[S{y,y)]), (14) 

D,ix,y) = (Tr[755(x,x)]Tr[755(t/,y)]), (15) 

where S{x,y) is the quark propagator. We calculate these disconnected contributions using 
the eigenmode decomposition of the quark propagator (ITT|) . Only the low-lying modes that 
we calculated for the Dirac operator spectral density are included in this analysis. We 
checked that the correlators at the long distances are unchanged when the number of low- 
modes are reduced from 50 to 40. Short-distance part is of course largely changed. 

The correlation functions of the relevant mesonic channels for am = 0.05, 0.025 and 0.01 
are plotted in Figure |6] from top to bottom panels. The correlators are for spatial directions 
and averaged over other directions. Left panels show the plots at /3 = 2.20, slightly below 
the transition temperature. We find clear distinctions among different channels at am = 
0.05 and 0.025. As the chiral limit is approached {am = 0.01), on the other hand, all 
the four channels are nearly degenerate. Since the disappearance of the near-zero modes 
becomes significant at am, = 0.01 (see Figure El middle panel), we infer that the difference 
among the four mesonic channels indeed originates from the near- zero eigenmodes of the 
Dirac operator. 

Above the transition temperature, /3 = 2.30 (T ~ 210 MeV), we find similar degeneracy 
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FIG. 6: Meson correlators at /3 = 2.20 (T ~ 180 MeV) and /3 = 2.30 (T ~ 208 MeV). Sea quark 
masses are am = 0.05, 0.025 and 0.01. Results for the vr, 6, rj and a channels are shown. Bands 
represent the statistical error. 
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FIG. 7: Same as Figure El but for /? = 2.40. 




FIG. 8: (Left panel) Contribution from the disconnected diagram to the iso-singlet pseudo-scalar 
correlator. Results at three different quark masses at /? = 2.30 are plotted. (Right panel) Relative 
contribution of the disconnected diagram to the connected is plotted in a logarithmic scale for some 
values of |x|. It shows that the disconnected contribution rapidly vanishes in the chiral limit. 

for small quark masses, as shown in right panels of Figure El At even higher temperature, 
/3 = 2.40 (T ~ 240 MeV), the degeneracy is found at higher quark masses (see Figure [7]). 
These observations are consistent with our interpretation that the near-zero eigenmodes 
below, say, 20 MeV are responsible for the splitting of the chiral partners for both iso-singlet 
and iso-triplet symmetries. 

Figure [H shows the contribution of the disconnected diagram at /3 = 2.30. It is clear from 
this plot that the disconnected contribution vanishes at smaller quark masses. 
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V. CONCLUSIONS 

Finite temperature phase transition of two-flavor QCD could be more complicated 
than previously thought due to the possibility of effective restoration of the axial U{1) 
symmetry. The standard analysis J3| assuming the pattern of symmetry breaking from 
SU{2)l ® SU{2)r (g) U{l)v to SU{2)v ® U{l)v can not be directly applied if the symmetry 
is effectively extended to include U{1)a- 

This work is one of the first attempts to understand the situation. Since the axial- 
anomaly sector of QCD is concerned, the standard technique of using the staggered fermions 
in lattice QCD is not appropriate for this purpose unless one approaches sufficiently close to 
the continuum limit. We instead use the overlap fermion formulation that exactly preserves 
chiral symmetry at the classical level and the U{1)a is violated by the axial anomaly as in 
the continuum theory. 

A side-effect of using the exactly chiral fermion formulation is that the global topology 
has to be fixed in the Monte Carlo simulation. Theoretical formula have been developed to 
obtain the correct physics of the 6 vacuum from such topology-fixed simulations. Essentially 
the error due to fixing topology is a finite volume effect of 0{1/V) and the leading 1/V 
correction can be extracted from the lattice data. We numerically tested the formula using 
quenched QCD at finite temperature, by calculating the topological susceptibility on the 
fixed-topology configurations and confirmed that the result agrees with that of the standard 
method. 

Having established that fixing topology at finite temperature does not introduce unex- 
pected systematic effects, we carried out a series of two-fiavor QCD simulations on 16^ x 8 
lattices around the transition temperature. By measuring the eigenmodes of the overlap- 
Dirac operator we found a gap of the spectral density in the chiral limit at high temperature, 
i.e. the near-zero modes disappear not just A ^ but including the modes of several tens of 
MeV. The disappearance of such near-zero modes has a correspondence with the degeneracy 
of meson correlators in the channels related by the U{1)a transformation, such as those be- 
tween TT and 5 oi 7] and a. The origin of this degeneracy is understood: the near- zero modes 
are the dominant source of the disconnected-diagram contributions to the meson correlators 
at long distances. 

For more quantitative understanding, several sources of systematic errors are to be ad- 
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dressed in the future works. Those include the finite volume effect and chiral extrapolation 
especially near the transition point. A cross-checking of the current results with an action 
that allows topology change while retaining a very good chiral symmetry would also be use- 
ful. Along this line, precise study of the order and critical exponents of finite temperature 
phase transition of two-fiavor QCD, which has not been established yet, will become feasible 
taking account of the effect of axial anomaly. 
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